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Abstract Using molecular dynamics computer simulations we investigate 
the out-of-equilibrium dynamics of a Lennard- Jones system after a quench 
from a high temperature to one below the glass transition temperature. By 
studying the radial distribution function we give evidence that during the 
^ • aging the system is very close to the critical surface of mode-coupling theory. 

Furthermore we show that two-time correlation functions show a strong de- 
'^ ', pendence on the waiting time since the quench and that their shape is very 

^ I different from the one in equilibrium. By investigating the temperature and 

Cj \ time dependence of the frequency distribution of the normal modes we show 

that the energy of the inherent structures can be used to define an effective 
(time dependent) temperature of the aging system. 

O 1 Introduction 
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G^ , In the last few years ample evidence has been accumulated that the mode- 

coupling theory of the glass transition (MCT) gives a reliable description of 
the dynamics of simple supercooled liquids on a qualitative as well as quanti- 
tative level [jl|. Recently it has even been documented that also some aspects 
of the dynamics of strong glass formers are described well by the theory [§]. 
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O ■ Thus we can conclude that many of the key aspects of the dynamics of su- 

^ ' percooled liquids are understood in a quite satisfactory way. This is not yet 

the case for the dynamics of glasses below the glass transition temperature, 
/\ • i.e. in that temperature regime in which the equilibrium relaxation time sig- 

j^ ■ nificantly exceeds the time scale of the experiment. Only relatively recently 

first attempts have been made to understand this out-of-equilibrium dynam- 
ics within the framework of statistical mechanics and thermodynamics [Q]. 
In particular it was found that, for certain systems, the equations of motion 
describing the dynamics below Tg are formally quite similar to the MCT 
equations, which, as discussed above, describe well the relaxation dynamics 
above Tg. Whether or not these out-of-equilibrium theories will be equally 
successful to describe the dynamics of structural glasses below Tg is currently 

"'^Talk presented at "Unifying Concepts in Glass Physics" , ICTP, Trieste 15 - 18 Septem- 
ber 1999 



not known and in the present paper we discuss some computer simulations 
which have been done to test these theories. 



2 Model and Details of the Simulations 

The system we consider is a binary (80:20) mixture of particles which interact 
with a Lennard- Jones potential, Vafs = '^^apii^'aplrY'^ — {(Tap/r)^]- Here 
a, P & {A, B} denote the species of the particles and the parameters ec/^ 
and cra/3 are given by eaa = 1-0, aAA = 1-0, bab = 1-5, (Tab = 0.8, esB = 
0.5, and a^B = 0.88. This potential is truncated and shifted at a distance 
aai3- In the following we will use cr aa and eAA as the unit of length and 
energy, respectively (setting the Boltzmann constant k-^ = 1.0). Time will 



be measured in units of Jma\p^/4:8eAA, where m is the mass of the particles. 

The equations of motions have been integrated with the velocity form 
of the Verlet algorithm, using a step size of 0.02. The number of A and B 
particles were 800 and 200, respectively and the size of the box was (9.4)^. 

In the past the equilibrium dynamics of this system has been determined 
in great detail [^ ^. In particular it was shown that at low temperatures the 
relaxation dynamics is described very well by MCT with a critical temper- 
ature of Tc = 0.435. To investigate the aging dynamics we therefore equili- 
brated the system at the high initial temperature Tj = 5.0 and quenched it at 
time zero to a final temperature Tf G {0.1,0.2,0.3,0.4,0.435}. This quench 
was done by coupling the particles every 50 time steps to a stochastic heat 
bath which was kept on during the subsequent propagation of the system at 
low (kinetic) temperature. In order to improve the statistics of the results 
we averaged over 8-10 different realizations of the system. 

3 Results 

Within the framework of the idealized version of MCT the aging process is 
viewed as a slow approach of the system to the so-called "critical surface" of 
MCT. This surface is a hyper-surface in the parameter space of the coupling 
constants, which in the case of a structural glass are given by the magnitude 
of S{q), the static structure factor at wave- vector q, and divides this space 
into a region in which the system is liquid-like and one in which it is solid-like, 
i.e. a glass. (In order to avoid some mathematical subtleties we consider only 
a discrete and finite set of wave-vectors, thus the parameter space is finite 
dimensional.) In order to check whether or not this surface does indeed 
have any relevance for the aging dynamics of our system we calculated the 
time dependence of gAA{f), the radial distribution function between two A 
particles, a quantity which is closely related to the static structure factor. In 
Fig. ^we show ^^aa for different times t after the quench (main figure). From 
this graph we see that immediately after the quench the function changes 
rapidly (compare the curves for t = and t = 10) but soon afterwards 



shows only a very weak time dependence and then can soon be considered as 
constant within the accuracy of the data. That this hmiting curve depends 
on the final temperature Tf of the quench is shown in the inset of the figure, 
where we show QAAif) for other values of Tf. We see that with decreasing 
Tf the height of the main peak increases and that its width decreases. The 
reason for this is that at such low temperatures the particles vibrate in the 
cages formed by their neighbors and that the size of this cage decreases 
with decreasing temperature. Note that this dependence demonstrates that 
for different values of Tf the system populates states in different regions in 
configuration space. 




Figure 1: Main figure: radial distribution function between the A particles 
for different times for Tf = 0.4 . The times are t = 0, i.e. before the quench, 
t = 10, 100, 1000, 10000, and 63100 time units. Inset: The same quantity at 
t = 63100 for different values of Tf. 



If the view put forward by MCT is correct these states should all be 
close to the critical surface discussed above. In order to check this pre- 
diction we calculated the area under the first peak in gAA{f). This area 
is roughly proportional to the height of the main peak in the static struc- 
ture factor and previous calculations have shown P, ^ that for simple sys- 
tems such as the present one this is the most relevant coupling parameter, 
i.e. the most relevant direction in the parameter space of the coupling con- 
stants. In Fig. ^we show this area c(t) for the different temperatures Tj, i.e. 
c(t) = An Jq" r'^gAA{t)dr, where Tc is the location of the minimum between 
the first and second peak in gAA^r). (Note that c is nothing else than the 
average number of A particles that surround an A particle.). The different 
symbols correspond to times 0, 10, 40, 60, 100, 160, 250, 400, 630, 1000, 
1580, 2510, 3980, 6310, 10000, 15850, 25120, 39810, and 63100 and are thus 



space roughly equidistant on a logarithmic time axis. (In order to expand the 
axis at low temperatures we chose a logarithmic temperature scale.) From 
the main figure we recognize that at the beginning the area does indeed de- 
pend strongly on time that, however, at intermediate and long times c{t) 
is essentially constant within the noise of the data (see inset). We see that 
within the time scale covered the value of c{t) changes only by about 6%. 
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Figure 2: Main figure: time dependence of the area under the first peak in 
5'aa('") for different values of Tf. The right-most curve is the same quantity 
for the equilibrium case. Inset: enlargement of the equilibrium curve around 
the critical temperature of MCT. 



That this small change is nevertheless very significant is demonstrated by 
the right-most curve (filled circles) in which the equilibrium value of c at 
different temperatures is shown [^ |TB[. From that curve we see that in the 
temperature range 5.0 >T> 0.446 the area changes also only by about 6%, 
despite the fact that the dynamics of the system slows down by about five 
orders of magnitude 0. 

The most relevant information from this graph is that the value of c{t) at 
long times seems to be almost independent of Tf, thus giving evidence that 
at these times the system is indeed close to the critical surface of MCT. That 
this is indeed the critical surface can be seen from the inset where we show 
an enlargement of the equilibrium curve at low temperatures. From previous 
investigations we know that at these temperatures the (equilibrium) system 
is very close to its critical temperature and that therefore the values of c of 
the equilibrium curve are very close to the critical ones and which are around 
11.93. From the inset we recognize that also the long time values of c{t) are 
very close to this number which is thus evidence that also the aging systems 
are, at long times, quite close to the critical surface. A careful inspection of 



the main figure reveals, however, that the curves for low values of Tf are a 
bit below this critical value, an observation which we will discuss below. 

From Figs. |l| and |^ we see that during the aging process the time de- 
pendence of the radial distribution function is rather weak. This situation 
is typical for so-called "one-time quantities", i.e. observables which in equi- 
librium are constant. (Below we will discuss exceptions to this trend.) A 
much stronger time dependence is found for the so-called "two-time quan- 
tities" , i.e. the generalizations of the equilibrium time correlation functions 
to the case of the out-of equilibrium case. In equilibrium a time auto- 
correlation function of an observable y{t) depends only on the time difference, 
i.e. {y{tu,)y{ty^ + r)) = {y{T)y{0)), where (.) is the thermodynamic average. 
This equahty no longer holds for the out-of equilibrium case, since due to the 
generation of the out-of equilibrium situation the time-translation invariance 
of the system is lost. Therefore it is necessary to keep track of both times, 
tui, the time since the quench, and r, the time since the start of the measure- 
ment. In the following we will study the case that the observable is ps{k,t), 
the space-Fourier transform of the density of a tagged at wave-vector k. This 
quantity is related to the positions of the particles via 

1 ^ 
^Pik,t) = -Y.exp[zk-rj{t)] . (1) 

^^ j=i 

In equilibrium the resulting time correlation function is the so-called inco- 
herent intermediate scattering function Fs{k, t) = {p{k, t)p{—k, 0)) which can 
be measured in scattering experiments. For the out-of equilibrium case we 
generalize this to 

Ck{tnj + ^,'tw) = j;tJ2 exp[«k ■ {vj{t^ + r) - Vj^r)] (2) 

j 

(Note that these last equations are trivially generalized to multi-component 
systems.) In Fig. |^ we show the r dependence of Ckit^ + T.t^) for the 
A particles for different waiting times t^. The value of the wave- vector is 
k = 7.23, the location of the main peak in the structure factor for the A-A 
correlation. (Other values of k, as well as the curves for the B particles look 
qualitatively similar.) From this figure we recognize that this time correlation 
function shows a very strong waiting time dependence, thus showing that 
the investigation of aging effects is much easier when one considers two-time 
quantities instead of one-time quantities (Figs. |l| and ||). Note that this 
waiting time dependence is not found at short times r, i/t^ is large, in that 
in this time regime the different curves collapse onto a master curve. In this 
time regime the particles are still inside the cage formed by their neighbors 
and thus we conclude from the figure that this vibrational motion becomes 
independent of the waiting time, if the latter is large. In Ref. IQ evidence was 
given that the approach of the curves to the plateau is given by a power-law, 
in agreement with the mean-field theories. 
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Figure 3: Time dependence of the generalization of the incoherent interme- 
diate scattering function to the out-of-equihbrium case. The different curves 
correspond to different waiting times. Tf = 0.4. The dashed curve is the 
equilibrium curve at 0.446. 



For larger times r the particle starts to leave the mentioned cage (the 
correlator starts to fall below the quasiplateau at intermediate times) and 
from the figure we see that the time at which this happens increases with 
increasing t^j in that the time scale for the second relaxation step increases 
rapidly with tw As it has been demonstrated elsewhere 0, this relaxation 
times scales roughly like t^^, i.e. it shows a sub-aging behavior. 

Also included in Fig. ^ is the equilibrium curve at T = 0.446 for the 
same wave- vector. Comparing this curve with the out-of-equilibrium curves 
for large t^ shows that the height of the plateau is very similar. What is, 
however, very different is the second relaxation process in that the equilibrium 
curve decays much more rapidly, i.e has a larger slope, than the aging curves. 
A detailed analysis shows that the former curve is approximated well by a 
Kohlrausch- Williams- Watts function ||^, exp(— (r/ro)^) whereas the latter 
curves are power-laws, with an exponent that depends on k but not on the 
waiting time [§, |l^ . 

The curves in Fig. ^ are for the final temperature Tf = 0.4, i.e. just 
about 10% below the critical temperature of MCT {Tc = 0.435). In order 
to see how Tf affects the relaxation we show in Fig. ^ the same type of 
correlation function as in Fig. ^ but this time for Tf = 0.1. Although the 
overall behavior of these curves are similar to the ones for the higher Tf some 
significant differences are found. First of all we see that the height of the 
plateau, fc{k), is now quite a bit larger than the one for Tf = 0.4. This 
change can easily be understood by recalling that for short times the motion 
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Figure 4: Time dependence of the generalization of the incoherent interme- 
diate scattering function to the out-of-equihbrium case for Tj = 0.1 



of the particles is dominated by their rattling inside the cage. To a first 
approximation this rattling can be described by a superposition of harmonic 
oscillators and thus their amplitude will be proportional to Tf. Thus we 
expect that 1 — fc{k) is proportional to Tf, and an inspection of the curves 
in Figs. ^ and ^ shows that this is indeed the case. More noteworthy is the 
observation that for large tw the curves seem to show a second plateau at 
long times, which seems not to be present for Tf = 0.1. 

In order to understand the origin of this second plateau it is useful to 
look at the individual runs for Tf = 0.1, which are shown for t,,. 
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the inset of Fig. ^. From this inset we see that the different curves show a 
relatively sharp drop in the time range 2 x 10^ < t < 2 x 10^ and then is 
almost constant. Note that the time at which this drop occurs depends on 
the realization. A careful analysis of the configurations just before and after 
a sudden drop shows that this fast relaxation is related to a very cooperative 
motion of about 10% of the particles [§. Thus it seems that this mechanism 
is the most effective way to release the stress that is in the configuration due 
to the quench. This is in contrast to the behavior at higher values of Tf in 
that there the stress is smaller and thus the system can remove it in a more 
gradual way, i.e. without the occurrence of the "earthquakes" that are seen 
at the lower temperatures. 

In the remaining part of the paper we will discuss the aging of the system 
from the point of view of the configuration space. For this we make use of 
the concept of the "inherent structure" (IS), which was introduced some time 



ago by Stillinger and Weber ||ll[ and can be described as follows: Any point 
in configuration space can be used as the starting point of a steepest descent 



procedure in the potential energy of the system. The endpoint of this steepest 
descent is the IS for the starting point. Thus in this way the configuration 
space can be decomposed uniquely into the basins of attractions of the IS 
(apart from some points of measure zero). By focusing on the IS we therefore 
can study the evolution of the system during the aging process without being 
disturbed by the vibrational part of the particle motion. 

In Fig. 1^ we show the temperature dependence of ejs the potential en- 
ergy of the system in the IS in equilibrium, Fig. |^a (see also Ref. |^). We 
recognize that at high temperatures eis is basically constant and starts to 
decrease quickly below T ^ 1.0 which shows that the energy landscape, as 
characterized by the height of the local minima, starts to change only when 
the system enters the supercooled regime. In Fig. Pd we show the time de- 
pendence of eis for the different final temperatures investigated. We see that 
the curves for small Tf show three regimes. (Although the curves with higher 
Tf show only two regimes we will argue below that also they should show 
a third regime if one would be able to simulate for longer times.) The first 
regime is observed at short times and in it eis(t) is essentially independent 
of time. After this time regime eis(t) enters the second regime during which 
the system is able to decreases it energy. After a certain time eis crosses over 
to a weaker time dependence and thus the system enters the third regime. 

In order to understand these two last time regimes it is useful to recall 
some results which have been obtained by analyzing the instantaneous nor- 
mal modes of supercooled liquids in equilibrium [|12| . Although these results 
have been obtained for supercooled water it is likely that they can be trans- 
ferred to the present systems as well. What has been shown in Ref. [|^ is 
that the number of modes that lead the system to a new local minimum de- 
creases with decreasing temperatures and vanishes at the MCT temperature 
Tc. Thus we can say that for temperatures above Tc the typical configura- 
tion of the system has at least one unstable mode, i.e. direction of motion, 
whereas for temperatures below T^ the system mainly sits in the vicinity of 
a local minimum, i.e. oscillates around a metastable location. Thus at T^ 
the system has a thermal energy which is comparable with the difference 
in energy between eis(Tc) and the lowest lying saddle point leading to a 
neighboring minimum. Hence, for temperatures below Tc the dynamics of 
the system becomes dominated by activated processes. Note that the fact 
that the temperature dependence of the relaxation dynamics shows a strong 
deviation from an Arrhenius temperature dependence |^, |n| implies that 
the height of the effective barrier between two minima is not constant but 
increases with decreasing temperature and, as explained above, at Tc this 
barrier is on the order of T^. 

All this holds for the equilibrium case. For the out-of-equilibrium case 
the situation is qualitatively similar but there is the important difference 
that now the system has only the thermal energy Tf. The three regimes 
seen in Fig. ^ can thus be explained as follows. In the first regime the 
typical configurations of the system are still quite similar to the ones found in 
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Figure 5: Time and time dependence of the energy of the IS. a): equihbrium, 
b) out-of-equihbrium. Reproduced from Ref. |TH . 



equilibrium at high temperatures, i.e. close to the initial temperature. In this 
part of configuration space the effective barriers between adjacent minima 
are relatively small but nevertheless noticeable. Since it takes the system 
some time to find configurations with lower energy, eis does not decrease. 
Only after some time the system manages to find configurations which are 
energetically more favorable and hence eis starts to decrease, i.e. in Fig. ^ 
eis(t) enters the second regime. The time until such better configurations 
are found increases with decreasing Ty, since a smaller kinetic energy makes 
it harder for the system to cross the barriers. 

During the aging process the system will lower its energy and start to 
explore configurations which, in equilibrium, correspond to lower and lower 
temperatures. After some time it will have reached that part of configura- 
tion space in which the effective barriers to cross from one minimum to the 
neighboring one have a height fcsT/- and thus the relaxation mechanism will 
change to an activated process. Since this type of relaxation is less efficient 
than the one in which the system still finds unstable modes, the rate by which 
the energy decreases is decreasing. Thus the eis(t) curve shows a bend, which 
can be seen in Fig. |^ when the system is entering the third regime. Note that 
this bend in the eis(t) curves should occur at a value of eis which increases 
with decreasing Tf and this is exactly what is seen in Fig. ^ In particular 
we expect from the reasoning above that if Tf is very close to Tc the third 
regime should hardly be visible, and this expectation is indeed supported by 
Fig. |. 

Since ejs seems to be a quite sensitive quantity to locate the place of the 
system in configuration space, and this is in contrast to most other one-time 



quantities, we can use it to define an effective temperature Te{t) during the 
aging process. For this we read off the value of eis(t) of the aging system at 
a time t, and define Te{t) to be that temperature T at which the system in 
equilibrium has the same value of eis (see Fig. |^). In order to check whether 
this definition of Te(t) has any physical meaning it is necessary to show that 
from the knowledge of T^it) it is possible to calculate other properties of the 
aging system. One such property is e.g. the distribution of the frequencies of 
the normal modes of the system. We have done this and found [|T^ that the 
value of Te{t) does indeed allow to calculate this distribution. In the present 
paper we show, however, only u, the first moment of this distribution (see 
Fig. I) since it has a better statistical accuracy than the distribution itself. 
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Figure 6: Time and time dependence of the first moment of the normal 
mode spectrum, a): equilibrium, b) out-of-equilibrium. Reproduced from 
Ref. 
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In Fig. P we show z/ for the equilibrium case as well as the aging case (left 
and right panel, respectively). A comparison of these curves with the one 
for eis in Fig. |^ shows immediately that the two sets of curves are very 
similar. Hence it follows that Te(t) can indeed be used to predict some of 
the properties of the aging system and thus can indeed be considered as an 
effective temperature. 

4 Summary and Conclusions 

We have investigated the properties of a simple glass former after a quench 
from a high temperature to a low temperature. From the radial distribution 
function we have evidence that for long times after the quench the system is 



in 



very close to the critical surface of the MCT. Note that, since this surface 
can be calculated from the structure factor, it will in the out- of- equilibrium 
situation be a function of Tf, since S{q) does not only depend on the IS but 
also on the thermal broadening of this configuration. This is the explanation 
why during the aging the system is "stuck" in parts of configuration space 
that depend on the value of Tf. At long times typical configurations have 
the property that if their IS is thermally broadened with a temperature 
Tf the resulting structure factor is very close to the critical surface of the 



MCT equations. Following the conclusions of Ref. |12] one can rephrase this 
by saying that for these configurations the typical barrier which leads to a 
neighboring minimum is on the order of Tf. 

Finally we mention that in a different place we have shown that the IS 
can be used to determine the configurational entropy of the system at low 



temperatures [|T6| which in turn allows to calculate such interesting quantities 
as the Kauzman temperature. 
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